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Abstract. 

We present detailed numerical calculations of the fast solar tachocline based on the assumption that the dynamo 
field dominates over the dynamics of the tachocline. In the present paper of the series, we focus on three shortfalls 
of the earlier models. First, instead of the simple oscillating dipole poloidal field we study the more general 
magnetic field structures reminiscen t of the butte r fly di agram. The migrating field is prescribed as the observed 
axisymmetric radial magnetic field (Stcnflo, Il988t ll994T). Our results are in good agreement with our analitical 
estimate and our previous works in lForgacs-Daika fc Petrovavl 11200 ll . I2002T) . but the polar "dip" in isorotational 
surfaces is strongly reduced in this case. On the other hand, a more realistic model should have a magnetic 
diffusivity decreasing significantly inside the radiative interior, so we also explore the effect of diffusivity and 
magnetic Prandtl number varying with depth. We found that the downwards decreasing magnetic diffusivity 
and Prandtl number have no significant effect on the solution, although the temporal variation of the tachocline 
thickness has decreased. 
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1. Introduction 

Helioseismic data allow us to study the structure and rota- 
tion rate of the solar interior. The tachocline is a thin tran- 
sitional layer below the solar convective zone. In this layer 
the surface-like rotation of the convection zone changes to 
near-uniform rotation in the radiative interior. The exis- 
tence and properties of this layer have be en known from 
helios eismic studies, recently reviewed bv ICorbard et alJ 
The tachocline is known to be extremely thin; 
the central radius of the tachocline appears to lie within 



to be < 0.05 1 


Kosovichevl Il99rl ICorbard et all 119991 


ISchou et all 11998 


: ICharbonneau et allll999l). 



Precise values for the parameters of the solar 
tachocline as such mean position, thickness, latitudinal 
and temporal variations, depend on the inversion tech- 
niques used for helioseismic data. lElliott fc Gouehl l)l999(l 
suggest that there exists a difference in the sound speed 
profile between the observations and the Standard Solar 
Models. By including an additional mixing layer below 
the convection zone, they calibrate the thickness of the 
tachocline and they find that the mean tachocline thick- 
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ness is 0.019 Rq with a formal standard error of about 5%. 
The observed light-element abundances at the surface also 
suggest that the mixing pr ocess in the solar tachoc line is 
confined to a shallow layer l|Brun et al l ll999M2002lh 

Making use of the most recent helioseismic results one 
can investigate the temporal variation and the lati tudinal 
structure of the tachocline. iBasu fc Antial l)200l|) found 
that the tachocline is prolate, with the difference between 
the tachocline position at 0° and that at 60° latitude being 
ab out 0.02 R^. This is i n agre e ment with resu l ts obt ained 
bv ICharbonneau et all l)l999lh IBasu fc Antial l)200l|) also 
report an increase in thickness of the tachocline with lat- 
itude, but this increase is less significant, though still 
at the 3ct level. Regarding the dynamical aspect, one 
would expect temporal variations asso ciated with the so- 
lar cycle to appear in the tachocline llHowe et all 12 000; 
IBasu fc Antial l200lt ICorbard et all l200l|) . iHowe et all 
(2000) report that a large-scale oscillation may be tak- 
ing place in this layer with a period about 1.3 yr. It is not 
clear whether this per iod is associated wit h solar-cycle- 
related variations, and lBasu fc Antial l)200l|) did not find 
any periodic or systematic changes in rotation rate in the 
tachocline region. This fact may be the consequence of the 
insufficiency of the helioseismic data. 
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The extreme thinness of the tachocline implies a 
strongly anisotr o pic tr ansport of angular momentum 
dSpiegel fc Zahrl Il992l) . Several different mechanisms 
have been proposed for this transport, but it is now widely 
believed that the magnetic field is instrumental in its ori- 
gin. It is known that an oscillating magnetic field can- 
not propagate far into the radiative zone. The extent of 
the penetration of this magnetic field is the skin depth 
Wcycj where r] is the magnetic diffusivity 
and w cyc is the frequency of the cycle. On the one hand, 
for 77 10 8 cm 2 /s the dynamo field cannot penetrate the 
tachocline, and we can expect the tachocline to be per- 
vaded by the internal primordial field. On the other hand, 
for T) <; 10 9 cm 2 /s the tachocline dynamics should be gov- 
erned by the dynamo field. As the associated diffusive 
and Lorentz timescales are also very different, these two 
cases basically correspond to the case of "slow" and "fast" 
tachocline, dis cussed in the literature (see esp. Table I in 
lOilmanLETOfll) . 

In recent years the case of a slow tachocline 
has been inv estigated extensively by a number 
of authors. iRiidiger fc Kitchatinovl l|l997|) and 
iMacGregor fc Charbonneaul l)l999j) studied the in- 
teraction between a large-scale field with fixed poloidal 
component and differential rotation without taking into 
consideration the meridional circulation. They found 
that the internal poloidal field of 10 -3 G is sufficient to 
confine the tachocline to its observed thickness if the 
internal field is fully conta ined within the radiative zone. 
iGough fc Mclntvref l|l998l) presented a model for the 
solar tachocline which allows the nonlinear interaction 
between the meridional flows and a la r ge-sca le magnetic 
field in self-consistent way. iGaraudl l|200l|) performed 
calculations taking into account the meridional flow and 
the self-consistent evolution of the poloidal field. 

In the first pa p er of these series 
ijForgacs-Daika fc Petrovavl l2002l hereafter Paper 
I), we studied the alternative case of a fast, turbu- 
lent tachocline with a turbulent diffusivity of about 
r\ = 10 10 cm 2 /s and we presented detailed numerical 
calculations allowing for the self-consistent evolution of 
the dipolar poloidal field. It was found that a sufficiently 
strong oscillatory poloidal field with dipolar latitude 
dependence at the tachocline-convective zone boundary 
is able to confine the tachocline to its observed thickness. 
This is in good agreement with o ur analitical estimate in 
IForgacs-Daika fc Petrovavl l|200l[) . 

In this paper, instead of a dipolar field we consider the 
more realistic case of a migrating field, studying its effect 
on the dynamics of the solar tachocline. After formulat- 
ing the mathematical problem in Sect. 2, results of our 
model containing a migrating magnetic field will be pre- 
sented in Sect. 3 We also explore the effect of the merid- 
ional flow and the impact of varying the diffusivity with 
radius. Finally, Sect. 4 summarizes the main results. 



2. The model 

In order to describe the differential rotation and the evo- 
lution of the large-scale magnetic field in the solar interior 
it is useful to write the equation of motion and the hydro- 
magnetic induction equation in a frame rotating with the 
fixed internal rotation rate. Thes e equatio ns can be found 
in the first paper of these series l|Paper J . Sect. 2.1). 

2.1. Flow Field 

The flo w field was treated in the same manner as in 
IPaper j nevertheless, so me addit ional remarks are taken 
here. As we mentioned in lPaperl Sect. 3.4, the total flow 
field can be written as the sum of the differential rotation 
and meridional circulation. While we compute the evolu- 
tion of the angular velocity with time, the meridional cir- 
culation is prescribed in our region, because in this case 
we focus on the effect of any meridional circulation in the 
tachocline, not on the complex problem of meridional cir- 
culation in the convective zone. Recent helioseismic inver- 
sions have shown that the outer convective envelope is per- 
vaded by a ~ 20 m s -1 poleward meridional flow. Such in- 
versions currently do not provide accurate meridional flow 
information at greater depth, thus we turn for guidance to 
the numerical simulations of therm ally driven turbulen t 



convection re c ently carri ed out bv | Mie sch eif^l 
iKuker fc Stixl l)200lj) and iKitchatinov fc Riidigei 



2000), 



1995) 



These simulations generate a mean meridional circulation 
and in these the meridional flow penetrates into the stable 
interior. In an attempt to incorporate at least qualitatively 
these observational and computati onal resu lts, we used a 
model similar to that described in lPanerll 

The strong subadiabatic stratification of the solar inte- 
rior sets an upper limit to the meri dional flow. An elemen- 
tary estimate discus sed inlPaper 3 yields v m ~ 10 cm s _1 . 
On the other hand iMiesch et al.l l)2000|) found that de- 
pending on assumed values of some of the input parame- 
ters this pattern is a persistent equatorward or poleward 
circulation in the overshoot region. Accordingly, the flow 
parameters used in our calculations were chosen to pro- 
duce simple smooth one- and two-celled flow patterns 
obeying the amplitude constraint as before mentioned, 
while reproducing the observed flow speed near the sur- 
face. For the density, we used the following expression: 



P(r) = C 



R, 



7 



(1) 



where the values of the parameters used are C — 1.2, 
m = 2. and 7 = 0.9. The flow patterns are shown in 
IPaper 1 Figs. 9 and 10. In these cases, the speed of the 
horizontal motion is ~ 5cm/s in the upper part of the 
radiative interior. 



2.2. Poloidal Field 

It is well known that sunspot activity occurs in the form of 
waves with periods ~ 11 yr, which propagate from middle 
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Fig. 1. Illustration of the evolution of the vector potential. 
Solid lines represent positive A and dashed lines denote 
negative A. 

latitudes towar ds the solar equator. In the e xtended but- 
terfly diagram l|Makarov fc Sivaramarl H98^) . one can see 
the poleward migration of the magnetic fields indicated 
by the migration pattern of a number of tracers such as 
qu iescent prominences or the coronal green line. 

IStenflol l|l988L Il994h analyzed 33 yr of synoptic obser- 
vations of the Sun's magnetic field carried out daily at the 
Mount Wilson and Kitt Peak observatories over all solar 
latitudes and longitudes. The observed longitudinal mag- 
netic field has been converted to a radial magnetic field 
assuming that the field direction is on average radial in the 
layers where the field is measured. This provides an op- 
portunity to analyse the global modes of the Sun. IStenflol 
dl988lll994h expanded the radial magnetic field in spher- 
ical harmonics and he found that the zonal magnetic field 
pattern can be represented as a superposition of N dis- 
crete modes with purely sinusoidal time variations with 
frequency uj cyc — 2n/22 yr -1 . Consequently, the evolu- 
tionary pattern of the axisymmetric radial magnetic field 
can be written in the following form: 

N 

B r (t,x) = ^ \ a l\ C0S [^cyc (t - U)] Pl(x), (2) 

with 



8=1 



<E>z = -uj cyc ti, 



(3) 



where B r is the radial, axisymmetric magnetic field, \ai\ 
is the amplitude, is the phase lag, U is the time lag and 

x = cos 9. 

IPetrovav fc Szakalvl (|l999j) found that the latitudinal 
distribution of the field at the surface reflects the con- 
ditions at the bottom of the convective zone, i.e. in this 
regard the convective zone behaves as a "steamy window" . 
Thus, if we assume that the above formula describes the 
time-evolution of the radial magnetic field at the base of 
the convective zone, then the development of the vector 
potential is given by: 



A{r,9,t) 



J B r (r, 9, t) sin 9rd9 + F{r, t) 



(4) 



where F is an arbitrary function. IStenflol l|l988l Il994j) also 
studied the dominance of the odd and even modes and he 



showed that the odd modes prevail in the evolution of the 
radial, axisymmetric magnetic field. 

Thus, the development of the vector potential is repre- 
sented as a superposition of the sinusoidal, 22 yr variations 
for the 7 odd modes with I = 1, 3, . . . , 13, as 



6 

A(9,t) = A ^2 V 1+2 k{9)\ai+2k\cos [uj cyc (t - t 1+2 k)} 

k=0 

Vt = lsin9 



2 

1 5 , 

V a = - sin 9 — - sin I 



7 21 

V 5 = ^ sin 9 - - sin 3 9 + — sin 5 ' 



V? — - sin ( 
2 



Vg 



1 . 
— sin 
2 
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4 

27 , 
IT 
11 



16 
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1 sin 
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(5) 
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256 
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256 
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■ sin 



8075 
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The factor Aq fixes the poloidal field amplitude. The \ai\ 
a mplitudes and_J>; phases of the odd modes were obtained 
bv IStenflol l|l994l) . The result is shown in Fig.Q] 



2.3. Equations for Our Model 

Based on |Papei"H Sect. 2.1, the azimuthal components of 
Eqs. (1-2), including the effects of diffusion, Coriolis force, 
meridional circulation, toroidal field production by shear, 
and Lorentz force, read 

d t uj = d r uj \ d r v + 4- + d r p- ) + decJ ^ C ° S n + (6) 
\ r p J sm 6 

QjW + flgw-^+L + C + Mi 



L = 



B 



(do A sin 9 — d r Ar cos ( 
(d e A&uv9 - Acos 9) 4 



47T/9r 3 sin 2 9 

d r B 
Airpr 2 sin 2 9 

9eB „ (—A sin 9 - d r Ar sin 9) 

3 sin 2 /) ' 

2^o cos 2 9 



Airpr 2, sin 2 9 



20 

C = ?A^(l-3cos 2 #)+d r ^ 

r J p r^p 

2 

Mi = u [ip (1 - 3 cos 2 9) + d r ijjr cos 2 9] + 
,1 — 3 cos 2 9 sin cos 9 

O r Wip 5 h OeOJOr'tp- 



r 2 p 



r 2p 
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Fig. 2. Spreading of the differential rotation into the radiative interior for rj = v = 10 10 cm 2 /s. Left-hand panel: 
contours of the time-average of the rotation rate under one dynamo period, with contour spacing 8.75 nHz; the solid 
lines represent the regions rotating faster and the dotted lines show the regions slower than £Io/2tt = 437nHz. Right- 
hand panel: the rotation rate at different latitudes as a function of radius. The peak amplitude of the poloidal magnetic 
field and the toroidal field are B p ~ 2400 G and B ~ 12000 G, respectively. 



,d r r) cos 2 6> 
o t B = B | j] 



r r 2 sin 2 I 



/ 2n „ \ „ n COS 9 
d r B — + d rV + d e B-! r — + 
\ r J r z sm 9 

d 2 r B V + 2 B^ + 

d r uj (do A sin 6 + A cos 9) + 

dgu) (-A sin 9 - d r Ar sin 9) -r + M 2 



M 2 = B 



?/)(3cos 2 6»- 1) 
1 



1 d r p 



d r il)cos 2 9] — 



d r B (1 - 3 cos 2 9) At- + 

sin cos 
d e B d r ip 

r Zp 



(7) 



* 0.04 
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GONG 




40 

Oolotitude 



Fig. 3. Latitudinal variation of the tachocline thickness 
for the case in Fig.|2 The solid line represents the thick- 
ness from the model and the dashed and dashed-dot 
lines show w from M DI and GONG data respectively 
l(Basu fc Antial l200ll) . The grey and the hatched areas 
show the error estimates on these. 



8tA 



M, = 



A— 

r 2 s 

d 2 A v 
A[il>(l- 
d r A^ (1 



d r A 



2 1) 



d„ A 



rj cos 9 



3 cos 2 



M 3 



3 cos 2 



d r ipr cos 2 9] -r — h 
r i p 

sin 9 cos 9 



1 

r 2 p 



(8) 



do A- 



r 2 p 



-d r ijj 



(For t he sake o f completeness we repeated here Eqs. (13 — 
15) of Paper I which contain a typo.) In these equations, 
L represents the terms associated with the Lorentz force, 
C denotes the terms associated with the Coriolis force 
and Mi, M 2 , M3 denote the terms associated with the 
advection by meridional circulation. 

2.4. Boundary Conditions for the Physical Parameters 

Our model does not include the convection zone, so the 
computational domain for the present calculations con- 
sists of just an annular meridional cut in the northern 




20 40 60 

Colotitude 9 [°] 

Fig. 4. Latitudinal variation of the radial jump in 57 for 
the case in Fig.|5| The solid line represents the thick- 
ness from the model and the dashed and dashed-dot 
lines show from MDI and GONG data respectively 
l|Basu fc Antial l200ll) . The dotted lines show the error es- 
timates on these. 
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hemisphere (0 < 9 < vr/2) in the upper part of the ra- 
diative interior, between radii r; n and rb cz , where rb cz is 
the radius of the bottom of the convection zone. For the 
solution of Eqs. (§HH) we use the following boundary con- 
ditions. 

Assuming the deep radiative core to be a rigidly ro- 
tating perfect conductor, we set both B and A to zero at 
the lower boundary. Both these quantities are also set to 
zero along the symmetry axis and they are antisymmetric 
about the equator. Furthermore, our upper boundary con- 
dition on the poloidal field is the solution of Eq. (0J. For 
the upper boundary condition chosen for the toroidal field 
we assume that the toroidal field is negligible in the con- 
vective zone compared to its value in the tachocline. Such 
a hypothetical situation is in line with current thinking 
in dynamo theory, and it may be a natural consequence 
of buoyancy-driven instabilities effectively removing any 
toroidal flux from the convective zone. 

For us we require axial and equatorial symmetry, and 
at the bottom of the convection zone we suppose that the 
rotation rate can be described with the same expression 
as in the upper part of the convection zone. Since the 
observations indicate that the differential rotation varies 
little within the convection zone and in the deep radia- 
tive interior, we can take d r us = at the bottom of the 
computational domain, so the lower boundary condition 
for us is different than in lForgacs-Daika fc Petrovavl l|200lL 
120021) where us — at the bottom of our box. Finally, we 
summarize and formulate these conditions in this way: 



4>= 0.00° (Equotor) 
0=15.71° 

= 60.00° 




t [years] 



Fig. 5. The thickness of the tachocline at different lati- 
tudes as a function of time for the case in Fig.|3 



The initial conditions chosen for all calculations are 



A(r,9,t = 0) = A J2 k=0 V 1+2k \a 1+2k \ cos($ 1+2fc ) 

at r — r^cz 

A(r, 0, t = 0) = at r < r bcz 

uj(r, 9,t = 0) = f^bcz - fio at r = r bcz 

Lj(r, 9, t = 0) = at r < rb cz 

B(r,9,t = 0) =0 



(11) 



Starting from the initial conditions, the solution is allowed 
to evolve in time until it relaxes to a very nearly periodic 
behaviour. 



A = A J2 k =Q V Wk\ai+2k\ cos [w cyc (t - ii +2 fc)] 

at r = rbcz 



US — fibc Z — 

B =0 

B = A = d r uj = 

B = A = d B uj = 

B =d e A= d e oJ = 



at r — ?"bc Z 
at r = rbc Z 
at r = r ln 
at 9 = 0, 
at 9 = tt/2, 



(9) 



where Qq is the rigid rotation speed at the observed rate 
Qq/2tt = 437 nHz and £! bcz is the rotation rate, which 
can be described with the same expression as in the upper 
part of the convection zone. 

In accordance with the observations of the GONG net- 
work, the expression used for S!b cz is 



2tt 



= 456 - 72 cos 2 9 - 42 cos 4 9 nHz. 



(10) 



2.5. Solution procedure 



The Navier-Stokes, the induction equations and the equa- 
tion for the poloidal magnetic field are coupled, nonsep- 
arable partial differential equations. We solve them as an 
initial-boundary value problem. We use a time relaxation 
method with a finite difference scheme first order accurate 
in time to solve the equations. A uniformly spaced grid is 
set up with 128 grid points in the r direction and 64 grid 
points in the 9 direction. 



3. Results 

As stated above, we examine the influence of a migrating 
field on the radial spreading of the differential rotation 
into the radiative interior. Prescribing a migrating mag- 
netic field at the top we define the form of the vector 
potential from the observed radial magnetic field. Besides 
the effect of the field strength, we study the consequences 
of the varying the diffusivity and the magnetic Prandtl 
number with radius, and the influence of the meridional 
circulation. 




I 



-4.00 u 

-8.00 

-12.00 



Time [yr] 



Fig. 6. Time-latitude diagram for the toroidal field 
(shaded region) and the radial field (contours) at the base 
of the convection zone. Equidistant contour levels of the 
radial field are separated by intervals of 30 G. 
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Fig. 7. Same as in Fig-121 but including the meridional circulation. In left-hand panel the circulation is counterclockwise, 
while in right-hand panel it is clockwise. 



3.1. Solution Without Meridional Circulation 

In the first calculations, only the field strength varies 
from case to case. The values of the viscosity and mag- 
netic diffusivity chosen for these simulations are identical, 
rj = v = I0 lo cm 2 /s, i.e. the magnetic Prandtl number 
is one Pr m = v/rj = I and the meridional circulation 
is neglected. Several calculations were run with different 
field strengths and we chose the most suitable case for us, 
where the sufficiently strong magnetic field is able to in- 
terfere the spreading of the differential rotation into the 
radiative interior and to reproduce the observed thickness 
of the tachocline. This case is shown in Fig.El after relax- 
ation 1 . 

In the left-hand panel we plotted the contours of the 
time-average of the rotation rate, which is defined as 

H(r, 6) = u(r, 9) + O (12) 

1 ft+Pcyc 

U(r,6) = — / u>(r,6,t)dt, (13) 

-'eye Jt 

where P cyc = 22yr is the period of the cycle. 

In the right-hand panel we plotted the rotation rate at 
di fferent la titudes as a function of radius. We note that 
in IPaper J we plotted the differential rotation amplitude 
Aw in right-hand panels, where u> was weighted with lat - 
itude and Aw was normalized (see Eq. (23) in IPaperll) . 
However, in this paper we retained the values of the fre- 
quency, so in this case the results can be compared easier 
with helioseismic observations and even the more insignif- 
icant variations of the rotation rate are more conspicouos. 
Based on the helioseismic observations we expect that 
the differential rotation changes to the rigid rotation in 
a thin transition region. It is well visible that the pre- 
scribed differential rotation at the base of the convection 
zone changes to the near-uniform rotation in a thin layer. 

We also study in detail the dynamics of the tachocline 
and compare to the observations. Accordingly, the thick- 
ness of the tachocline and the change in the rotat ion rate 
across the tachocline are defined as in iBasul (|l99jj). Thus, 
the jump in the rotation rate across the tachocline 5tt is 
the difference between the rotation rate at the base of the 

1 Color figures and computer animations illustrat- 
ing the time development of some of our solutions 
can be downloaded from the following web site: 
http : / /astro . elte .hu/kut at /sol/fast 2/fast2e .html 



0.08 








0.06 


No circulation 

Equatorward ci 

Poleward circu 


culation 
zition 




0.04 




1 1 




0.02 


-Ju 
- \ \ I 
M 1 


'i 1 
1 1 u 

1 ;'V. 




0.00 








c 


20 40 

Colatitude 


60 

e [°] 


80 



Fig. 8. Latitudinal variation of the tachocline thickness. 
The solid line represents the case shown in Fig. El and the 
dashed and dashed-dot lines respectively denote the thick- 
ness including the meridional circulation. 

convection zone and the rotation rate in the interior. The 
thickness of the tachocline w is defined as the rotation rate 
increases from the factor 1/(1 + e) of its maximum value 
to the factor 1 — 1/(1 + e) of its maximum value in the 
range r = r t — w to r — r t + w. 

The fast tachocline shows a quite marked latitude- 
dependence in our results, as illustrated in Figs. 13141 In 
both figures we also plot the results o f the helioseis- 
mic m easurements from MDI and GONG l|Basu fc Antial . 
l200lh . It should be noted that, according to the helio- 
seismic results, the latitudinal variation in position of 
tachocline is significant, while the variation in width is 
not clear, but the results tend to suggest that the thick- 
ness increases with latitude. Since the top of the compu- 
tational domain is fixed at 0.69-R© in our model, therefore 
the position of the tachocline cannot be directly compared 
with the observations. Notice that the thickness of the 
tachocline is close to zero at around the colatitude of 60°, 
because at this colatitude the rotation rate is equal to the 
rigid rotation velocity. This region is the so-called corota- 
tion zone (Fig.^J). In addition to latitudinal variations, we 
also investigate the temporal variations in the properties 
of the tachocline, and we find that the thickness of the 
tachocline depends on cycle phase (Fig.[SJ). This variation 
is not confirmed by the observations, but the helioseis- 
mic data set is relatively short. The development of the 
toroidal magnetic field and the radial magnetic field shown 
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in Fig. correspond to the observed evolutionary pattern 
of the magnetic field, i.e. the butterfly diagram. 



layer to its observed thickness. We studied it for the sake 
of completeness. 



3.2. Solution With Weak Meridional Circulation 

In the second case, we study the influence of the merid- 
ional circulation. As we mentioned, the strong subadia- 
batic stratification sets an upper limit to the amplitude of 
the meridional circulation (v m ~ 10 cm s _1 ). 

We use the flow patterns shown in lPaper 1 Figs. 9 and 
10, where the speed of the horizontal flow is weak. In these 
cases, we considered the effect of the meridional flow on 
the magnetic field (M 2 , M 3 in Eqs. GEJl)- The results with 
the meridional circulation included can be seen in Fig.0 
The poloidal field amplitude used here was the same as 
in the previous case. Notice that the effect of the weak 
meridional flow of ~ 5 cm s _1 on the magnetic field is not 
dramatic. 

In Fig.[7] it can be seen that the rigid rotation 
rate is shifted compared to the first case in accor- 
dance with the direction of flow. An equatorward flow in 
the tachocline region increases the near-uniform rotation 
rate, while a poleward flow reduces it. We remark that 
the lower boundary condition fo r ui is different than in 
iForgacs-Daika fc Petrovavl (|200lL 12002^ where u = at 
the bottom of our box. In these cases d r uj = at the 
bottom of the computational domain, so the effect of the 
meridional circulation on the rotation is perceptible in the 
radiative interior. It is apparent that the deep meridional 
flow changes the rigid rotation rate. If we raised the am- 
plitude of the meridional circulation at the upper part 
of the radiative interior (v m ~ 2 m s _1 ), but we reduced 
the penetration of the meridional circulation the merid- 
ional motion is in the overshoot region -, the shift is not 
dramatic. This result could play an important role in the 
model for the solar dynamo models wherein a meridional 
flow penetrates into the stab le layers below the tachocline 
iNandvfc ChoudhurilEooi. 

In Fig.[S|the differences between the three results are 
apparent. Note that the rotation rate in the interior is 
changed in these cases and we used it to define the thick- 
ness of the tachocline. It is well visible that the effect of 
the weak circulation is not dramatic and chiefly limited to 
the pole. In this paper the toroidal field strength is higher 
as in the earlier paper and on the other hand the upper 
boundary condition on the poloidal field is controlling the 
influence of the meridional motion. In case of equatorward 
circulation it can be seen that the tachocline thickness at 
the pole extended and in the corotation zone the peak is 
apparent. Both are due to the transport by the equator- 
ward motion. If we assume that in the tachocline region 
the circulation is toward the pole then the thickness at 
the high latitudes is reducible. It must be noted of course 
that the solar dynamo theories favour the counterclock- 
wise circulation. On the other hand in our model we have 
not emphasized the role of the circulation, but the role of 
the Maxwell stresses which are capable to limit the shear 



3.3. Varying the diffusivity with radius and Pr m = 1.0 

In this section t he diffusivity is allowed to vary with ra- 
dius. In lPaperll we only investigated cases where the dif- 
fusivities are constant throughout the computational do- 
main. The question may arise that the high viscosity in 
the deep radiative zone could be a disfigurement to our 
results? In order to answer, we change the magnetic dif- 
fusivity and the viscosity with radius as 



r) (or v) = r/o (or vq) 



1 + exp r 



fbcz 



(14) 



where 770, vo is the maximum value of the diffusivity and 
viscosity, T and r v are geometric parameters. Important 
notice that the used numerical technique sets a limit to 
this variation with radius. Namely, the lowest value of vis- 
cosity or magnetic difusivity in the domain determines the 
diffusive timescale as against the time step which is deter- 
mined by the highest value of these. Hence the very low 
value of viscosity gives rise to the extremly long integra- 
tion times. 

Thus, in this case we set these parameters as follows: 
i]a (or u ) = 10 10 , T = 16 and r v = 4.2 x 10 10 , and the 
magnetic Prandtl number is one Pr m = v/rj = 1. The 
magnetic diffusivity is varied with depth between approx. 
10 5 and 10 10 . This variation of the diffusivity, while quan- 
titatively much milder than expected physically, provides 
an opportunity for us to consider the effects of a variable 
diffusivity without the necessity of impossibly long inte- 
gration times. 

In left-hand panel of Fig.Elit is clearly visible that the 
influence of the downwards decreasing diffusivity on the 
solution is small. However it is noticeable that in the po- 
lar region the differential rotation penetrates into deeper 
layers of the radiative interior. On the other hand the am- 
plitude of the temporal variation of the thickness is smaller 
than in the case of Fig.0 

3.4. Varying the magnetic Prandtl number with depth 

Under realistic conditions, the turbulent magnetic diffu- 
sivity in the tachocline should be reduced to the micro- 
scopic magnetic diffusivity in the radiative core and the 
turbulent viscosity to molecular and radiative viscosity. 
However the magnetic diffusivity is greater than the vis- 
cosity a nd it is less reduced than the viscosity (see esp. 
Fig. 1 in lRiidiger fc Kitchatinovl Il997|) . so the magnetic 
Prandtl number is varied with depth. 

As argued earlier, the low diffusivities often make the 
numerical simulations problematic, therefore we tried to 
give a qualitative test instead of a quantitative analysis. 
Accordingly, we used the above formula l|14f> to define the 
variation of rj and v with radius and the required param- 
eters are following: rj Q = 1 x 10 10 , T = 8, r v = 4.2 x 10 10 
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Fig. 9. Same as in Fig-EI Left-hand panel: the diffusivity is varied with radius. Right-hand panel: the magnetic Prandtl 
number is varied with depth. 




t [years] 



Fig. 10. The thickness of the tachocline at different lat- 
itudes as a function of time for the case in right-hand 
panel of Fig. (varying the magnetic Prandtl number with 
depth). 



for j] and u = 10 9 , T = 10, r v = 4.2 x 10 10 for v, so the 
magnetic Prandtl number is varied between 0.024 and 0.1. 

Our result is illustrated in right-hand panel of Fig.|U| In 
this case the poloidal field strength is B p ~ 1500 G and the 
maximum value of the toroidal field is B ~ 6000 G. It is 
well visible that the penetration of the differential rotation 
into the radiative interior is strongly reduced in the polar 
region and the temporal variation of the thickness is also 
decreased (see Fig lTT))) . 

4. Discussion and Conclusion 

The models discussed in this paper show that the dy- 
namo field dominates over the dynamics of the fast solar 
tachocline provided that the turbulent diffusivity in the 

tachocline region excee ds n <; 10 9 cm 2 /s. 

In our first work (|Forgacs-Daika fe Petrovavl l200l[) 
we presented an analytical estimate for the strength of 
an oscillatory poloidal field, which is able to confine the 
tachocline to its observed thickness for r) ~ 10 10 cm 2 /s. 
Based on the estimate we showed a simple numerical 
calculation. A crucial assumption in this calculation is 
that the oscillatory poloidal magnetic field does penetrate 
below the convective zone to depths comparable to the 
tachocline thickness and the penetration is prescribed a 
priori. 



In our forthcoming model ijPaper J) , instead of the a 
priori prescribed oscillatory poloidal magnetic field we in- 
cluded the evolution equation for the poloidal field, but 
the poloidal field was a simple oscillating dipole, which 
was not very realistic. In that work we also explored the 
parameter space a bit more widely by varying the mag- 
netic Prandtl number and the diffusivity, but we did not 
investigate the effects of the variations of these parame- 
ters with depth. We found that the confining field strength 
increased with increasing viscosity and magnetic Prandtl 
number. 

In the present work we focus on three shortfalls of 
the earlier models. First, instead of the simple oscillat- 
ing dipole poloidal field we study the more general mag- 
netic field structures reminiscent of the butterfly diagram. 
On the other hand, a more realistic model should have a 
magnetic diffusivity decreasing significantly inside the ra- 
diative interior. 

For the more general magneti c field structures, we use 
the results of lStenflol lll988ill994h . where the evolutionary 
pattern of the axisymmctric magnetic field can be written 
as a superposition of the sinusoidal, 22 yr variations for 
the 7 odd modes spherical harmonic components. Based 
on our previous works we studied the field strength and we 
found that a sufficiently strong migrating magnetic field is 
able to reproduce the observed thic kness of the tachocline. 
This agrees with the results in lForgacs-Daika fe Petrovavl 
apart from the fact that the time-latitude 
diagram for the toroidal field (Fig.EJ) is now similar to the 
butterfly diagra m. The lower boundary condition f or u> i s 
different than in lForgacs-Daika fc Petrovavl l|200ll 12002^ . 
giving us a possibility to study the effect of the dynamics 
of the upper layers on the inside layers. In Fig.[3it can be 
seen that because of the high diffusivity there is a little 
deviation in the rotation rate compared to the fixed in- 
ternal rotation rate Slo at different latitudes. Furthermore 
an equatorward meridional circulation in the tachocline 
region increases the near-uniform rotation rate, while a 
poleward flow reduces it. Note that this phenomenon is 
sensitive to the penetration depth of the flow into the up- 
per radiative interior. We found that the internal rotation 
rate is not shifted extremly in spite of the increasing am- 
plitude of circulation if the circulation is dominated in the 
upper part of the tachocline. The downwards decreasing 
diffusivity or rather the magnetic Prandtl number does 
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not appreciably affect the solution. Evidently, this model 
is simplified treatment of the turbulence in the tachocline, 
nevertheless it gave us an opportunity to consider these ef- 
fects. The varying magnetic diffusivity and Prandtl num- 
ber with depth give rise to a change of the temporal vari- 
ation of the tachocline thickness. It is shown in Fig.llUI 
that the downwards decreasing magnetic Prandtl number 
tends to reduce the spatiotemporal variation. 

The main results of present and earlier models is that 
the dynamo field in the turbulent tachocline by Maxwell 
stresses is capable to limit the shear layer to its observed 
thickness. The next step of this analysis consists in com- 
paring t he numer ical results in the present work and the 
model in IPaperll with the helioseismic results. In |Paper"H 
the thickness of the tachocline showed the "pit" at the 
pole, which is not confirmed by the available helioseismic 
information, which may be in conseque nce of the latitu- 
dinal limit of observations. In IPaper J the poloidal field 
was a simple oscillating dipole, and the reality of this pit 
needs to be verified with other poloidal field geometries 
and other parameters. In this paper we found that in case 
of the more general magnetic field structure the polar pit 
is reduced to its half and in case of varying the magnetic 
Prandtl number with depth the pit is almost disappeared 
fFig. lll(l . If we use the more general field geometry remi- 
niscent of the butterfly diagram, the penetration depth of 
the differential rotation into the radiative interior in the 
polar region is decreased owing to the fact that at higher 
latitudes the stonger magnetic fields stay up. In case of the 
decreasing magnetic Prandtl number the influence of the 
Maxwell stresses are in the ascendant compared with the 
viscous stresses, so the permeation of the differential rota- 
tion below the convective zone at the pole reduces further. 
On the other hand regarding the dynamical aspect, we 
studied the temporal variations of tachocline. Due to the 
relatively short helioseismic data set, the temporal vari- 
ations associated with the solar cycle to appear in the 
tachocline is not clear. In both models, the vector poten- 
tial of the magnetic field is prescribed with a period of 
22 years at the top of the computational domain, so we 
found that the thickness of the tachocline depends on cycle 
phase. However the amplitude of the temporal variations 
can be reduced if the impact of Maxwell stresses is en- 
hanced compared with the viscous stresses - for example 
the magnetic Prandtl number is less than one (cf . Fig. llQI) . 

As stated earlier, an important shortfall of these mod- 
els is obviously the extremely simplified treatment of the 
turbulent transport. If the turbulence arises from the in- 
stabilities in the tachocline, then the diffusivities should 
clearly be some, possibly very complex, functionals of the 
velocity field and of the magnetic field. The next step to 
model this would be to use a Smagorinsky-type formula 
or k — e model for the diffusivities. 

Other important points what we should consider in 
the future are the boundary conditions on the mag- 
netic field and the complete dynamo-generated poloidal 
field, because the conditions in the present treatment 
poorly represent the dynamical interaction between the 




Fig. 11. Latitudinal variation of the tachocline t hickness. 
The solid line represents the case shown in Fig. 6 of |Paper"3 
and the dashed and dashed-dot lines respectively denote 
the thickness including the meridional circulation. 

flow and the m agnetic field in the tachocline and the 
convection zone l|Petrovavl l200Cfl . Many theoretical mod- 
els suggest that the toroidal magnetic field is up to 
10 5 G in the tac hocline, required to produc e suns pots 
in low latitudes IjPetrovav fc Moreno-Insertist [1997) . In 
our model the peak toroidal field is above about 10 4 
G, but this may change if different boundary condi- 
tions are used or the dynamo-generated field is incl uded. 
iDiknati fc Oilman! ifeOOlh : lOilman fc Dikoatil l|2002l) sug- 
gest that the dynamo-generated toroidal field should play 
an important role in the global instability of tachocline 
differential rotation. In future work we plan to develop 
our model along the lines described above. 
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